!******************************************************************************************************************
!Module containing:
!(a) SUBROUTINE param_initialization: the initialization of the the parameter vector (global variable)
!(b) SUBROUTINE matrix_initial: the subroutine to initialize the matrix of the parameters for the amoeba routine   
!******************************************************************************************************************
MODULE mod_param_initial

USE mod_types


IMPLICIT NONE


!****************************
!(1) Declare variables
!****************************
INTEGER, PARAMETER:: tot_param = 143
! p_vector: it contains the vector of guesses
REAL(dp), DIMENSION(tot_param):: p_vector 

INTEGER, DIMENSION(tot_param):: est_ind
INTEGER, ALLOCATABLE:: pos(:)

INTEGER:: num_param

CONTAINS

!*******************************************************************************
!(2) Initialize the vector of parameters and specify which are to be estimated
!*******************************************************************************
SUBROUTINE param_initialization

IMPLICIT NONE


!outcome probabilities
REAL(dp):: alpha_1_e, beta_1_e, alpha_2_e, beta_2_e, alpha_3_e, beta_3_e, &
           alpha_11_e, beta_11_e, alpha_21_e, beta_21_e, alpha_31_e, beta_31_e, &
		   alpha_12_e, beta_12_e, alpha_22_e, beta_22_e, alpha_32_e, beta_32_e, &
		   alpha_13_e, beta_13_e, alpha_23_e, beta_23_e, alpha_33_e, beta_33_e, &
		   alpha_14_e, beta_14_e, alpha_24_e, beta_24_e, alpha_34_e, beta_34_e
! separation rate
REAL(dp):: xi_1_e, xi_2_e, xi_21_e, xi_3_e, xi_11_e,  xi_14_e, xi_24_e, xi_12_e, xi_32_e, &
           xi_33_e, xi_34_e, xi_15_e, xi_26_e
!outcome realizations
REAL(dp):: y_1_h_e, y_1_l_e, y_2_h_e, y_2_l_e, y_3_h_e, y_3_l_e, y_high_e, tau_e, tau, &
           output_edu_e, output_year_e, basel01_e, basel11_e,  basel21_e,  basel31_e, &
		   basel02_e, basel12_e,  basel22_e,  basel32_e, &
		   basel00_e, basel10_e,  basel20_e,  basel30_e, &
		   c01_e, c21_e, c31_e, c02_e, c12_e, c32_e, c03_e, c13_e, c23_e, cpen_e


REAL(dp):: B_21_e, B_22_e, B_23_e, B_11_e, B_12_e, B_13_e
REAL(dp):: PBAR_e, delta_e
REAL(dp):: gamma10_e, gamma11_e, gamma12_e, gamma13_e, gamma20_e, gamma21_e, gamma22_e, &
           gamma23_e, gamma30_e, gamma31_e, gamma32_e, gamma33_e, &
		   gammapast11_0_e, gammapast11_1_e, gammapast11_2_e, &
		   gammapast21_0_e, gammapast21_1_e, gammapast21_2_e, &
		   gammapast22_0_e, gammapast22_1_e, gammapast22_2_e, &
		   gammapast32_0_e, gammapast32_1_e, gammapast32_2_e, &
		   gammapast33_0_e, gammapast33_1_e, gammapast33_2_e, &
		   gammapast11_3_e, gammapast22_3_e, gammapast33_3_e, gammapast21_3_e, gammapast32_3_e, &
		   gammafut12_e, gammafut23_e

REAL(dp):: p_1_e, p_2_e, gammap_01_e, gammap_11_e, gammap_21_e, gammap_31_e,&
           gammap_02_e, gammap_12_e, gammap_22_e, gammap_32_e,&
		   gammap_03_e, gammap_13_e, gammap_23_e, gammap_33_e, gammap_41_e, gammap_51_e, &
		   gammap_42_e, gammap_52_e, gammap_43_e, gammap_53_e, p_3_e, gammap_04_e, &
		   gammap_14_e, gammap_24_e

REAL(dp):: p_70_e, p_71_e, p_72_e, p_73_e, p_74_e, p_75_e, p_76_e, p_77_e, p_78_e, p_79_e

REAL(dp):: s_0_e, s_1_e, s_2_e, s_3_e, sigma1_u_e, sigma12_u_e, sigma13_u_e, sigma14_u_e, sigma2_u_e, sigma3_u_e, sigma32_u_e, sigma33_u_e, sigma34_u_e, sigma_eps_e, sigma_vu_e, &
           price_e, s1_0_e, s1_1_e, s1_2_e, s1_3_e, sigma22_u_e, sigma23_u_e, sigma24_u_e, &
		   s_4_e, s1_4_e, s_5_e, s1_5_e, s_6_e, p1_low1_e, p1_low2_e, p1_low3_e, &
		   p1_high1_e, p1_high2_e, p1_high3_e , s_02_e, s_03_e, s_002_e, s_0022_e, s_0023_e, s_0024_e, &
		   s_003_e, s_0032_e, s_0033_e, s_0034_e, &
		   omega_0_e, omega_11_e, omega_21_e, &
		   omega_12_e, omega_22_e, omega_13_e, omega_23_e, omega_3_e, omega_4_e, lambda_11_e, lambda_12_e, lambda_13_e, &
		   s_12_e, s_22_e, s_13_e, s_23_e, a1_e, b1_e, a2_e, b2_e, a3_e, b3_e, tau2_e,tau3_e, &
           s_52_e, s_53_e, s_54_e, s1_52_e, s1_53_e, s1_54_e, d_e, d1_e, d2_e, d3_e, d4_e, d5_e, d6_e, d7_e, d8_e, d9_e, &
		   d10_e,  d11_e, d12_e, phk_e, phk2_e, phk3_e, d13_e, d14_e, d15_e, d16_e, d17_e, &
		   d1_1_e, d1_2_e, d1_3_e, xi_25_e, xi_35_e, d18_e, d19_e, d20_e, d21_e, d22_e, d23_e, d24_e, d25_e, d26_e, d27_e, &
		   d28_e, d52_e, d29_e, d30_e, d31_e, d32_e, plus_e, d5j2_e, d5j3_e, d53_e, d533_e, d523_e, wd_e, wd3_e, wd2_e, tencoef_e, &
		   tencoef1_e, tencoef2_e, wd52_e, d5t2_e, d5t3_e, d5t4_e, wpen_e, wvar_e, v74_e, v75_e, v76_e, v77_e, v78_e, v79_e, wpen1_e, &
		   por_e, prop_e, tencoef3_e, tencoef31_e, d5j22_e, d5j23_e, d5j24_e, tencoef22_e
 
REAL(dp):: w_70_e, w_71_e, w_72_e, w_73_e, w_74_e, w_75_e, w_76_e, w_77_e, w_78_e, w_79_e, w3_74_e, w3_75_e, w3_76_e, w3_77_e, w3_78_e, w3_79_e, &
           
		   w2_74_e, w2_75_e, w2_76_e, w2_77_e, w2_78_e, w2_79_e

REAL(dp):: s5_1_e, s5_2_e, s5_3_e, s1_55_e, d55_e, s56_1_e, s56_2_e, s56_3_e, s3_1_e, s3_2_e, s3_3_e


REAL(dp), PARAMETER:: scaleb = 100000.00_dp, scalea = 100000.00_dp


!--------------------------------
!(a) Read the parameter values 
!--------------------------------
OPEN(unit=8000, file = 'param_values.txt', position='rewind')

READ(8000,*) alpha_1_e
READ(8000,*) beta_1_e
READ(8000,*) alpha_2_e
READ(8000,*) beta_2_e
READ(8000,*) alpha_3_e
READ(8000,*) beta_3_e
READ(8000,*) xi_1_e
READ(8000,*) xi_2_e
READ(8000,*) xi_3_e
READ(8000,*) p_1_e
READ(8000,*) gammap_01_e
READ(8000,*) gammap_11_e
READ(8000,*) gammap_02_e
READ(8000,*) gammapast11_0_e
READ(8000,*) gammap_03_e
READ(8000,*) s_0_e
READ(8000,*) s_1_e
READ(8000,*) s_2_e 
READ(8000,*) sigma1_u_e
READ(8000,*) sigma_eps_e
READ(8000,*) s1_5_e
READ(8000,*) s_3_e
READ(8000,*) s_002_e
READ(8000,*) s_003_e
READ(8000,*) sigma2_u_e
READ(8000,*) sigma3_u_e
READ(8000,*) w_74_e
READ(8000,*) w_75_e
READ(8000,*) w_76_e
READ(8000,*) w_77_e
READ(8000,*) w_78_e
READ(8000,*) w_79_e
READ(8000,*) d2_e
READ(8000,*) d4_e
READ(8000,*) gammap_04_e
READ(8000,*) s1_52_e
READ(8000,*) s1_53_e
READ(8000,*) s1_54_e
READ(8000,*) d6_e
READ(8000,*) d9_e
READ(8000,*) d10_e
READ(8000,*) d11_e
READ(8000,*) xi_24_e
READ(8000,*) sigma32_u_e
READ(8000,*) sigma33_u_e
READ(8000,*) sigma34_u_e
READ(8000,*) d7_e
READ(8000,*) d8_e
READ(8000,*) d5_e
READ(8000,*) d1_e
READ(8000,*) xi_26_e
READ(8000,*) phk2_e
READ(8000,*) phk3_e
READ(8000,*) s_4_e
READ(8000,*) d12_e
READ(8000,*) xi_32_e
READ(8000,*) xi_11_e
READ(8000,*) xi_14_e
READ(8000,*) p_2_e
READ(8000,*) p_3_e
READ(8000,*) xi_34_e
READ(8000,*) xi_15_e
READ(8000,*) d13_e
READ(8000,*) d14_e
READ(8000,*) d15_e
READ(8000,*) d_e
READ(8000,*) d16_e
READ(8000,*) d17_e
READ(8000,*) xi_21_e
READ(8000,*) d1_1_e
READ(8000,*) d1_2_e
READ(8000,*) d1_3_e
READ(8000,*) y_1_h_e
READ(8000,*) cpen_e
READ(8000,*) gammapast11_1_e
READ(8000,*) xi_25_e
READ(8000,*) xi_33_e
READ(8000,*) sigma22_u_e
READ(8000,*) sigma23_u_e
READ(8000,*) sigma24_u_e
READ(8000,*) s_0022_e
READ(8000,*) s_0023_e
READ(8000,*) s_0024_e
READ(8000,*) s_0032_e
READ(8000,*) s_0033_e
READ(8000,*) s_0034_e
READ(8000,*) sigma12_u_e
READ(8000,*) sigma13_u_e
READ(8000,*) sigma14_u_e
READ(8000,*) gammapast33_0_e
READ(8000,*) gammapast33_1_e
READ(8000,*) d2_e
READ(8000,*) d18_e
READ(8000,*) d19_e
READ(8000,*) d20_e
READ(8000,*) d21_e
READ(8000,*) d22_e
READ(8000,*) d23_e
READ(8000,*) d24_e
READ(8000,*) d25_e
READ(8000,*) d26_e
READ(8000,*) d27_e
READ(8000,*) d28_e
READ(8000,*) d52_e
READ(8000,*) d29_e
READ(8000,*) d30_e
READ(8000,*) d31_e
READ(8000,*) delta_e
READ(8000,*) d_e
READ(8000,*) d32_e
READ(8000,*) plus_e
READ(8000,*) d5j2_e
READ(8000,*) d5j3_e
READ(8000,*) d53_e
READ(8000,*) d533_e
READ(8000,*) d523_e
READ(8000,*) wd_e
READ(8000,*) wd3_e
READ(8000,*) wd2_e
READ(8000,*) tencoef_e
READ(8000,*) tencoef1_e
READ(8000,*) tencoef2_e
READ(8000,*) wd52_e
READ(8000,*) d5t2_e
READ(8000,*) d5t3_e
READ(8000,*) d5t4_e
READ(8000,*) wpen_e
READ(8000,*) wvar_e
READ(8000,*) v74_e
READ(8000,*) v75_e
READ(8000,*) v76_e
READ(8000,*) v77_e
READ(8000,*) v78_e
READ(8000,*) v79_e
READ(8000,*) wpen1_e
READ(8000,*) por_e
READ(8000,*) prop_e
READ(8000,*) tencoef3_e
READ(8000,*) tencoef31_e
READ(8000,*) d5j22_e
READ(8000,*) d5j23_e
READ(8000,*) d5j24_e
READ(8000,*) w2_74_e
READ(8000,*) w2_75_e
READ(8000,*) w2_76_e
READ(8000,*) w2_77_e
READ(8000,*) w2_78_e
READ(8000,*) w2_79_e
READ(8000,*) w3_74_e
READ(8000,*) w3_75_e
READ(8000,*) w3_76_e
READ(8000,*) w3_77_e
READ(8000,*) w3_78_e
READ(8000,*) w3_79_e
READ(8000,*) s5_1_e
READ(8000,*) s5_2_e
READ(8000,*) s1_55_e
READ(8000,*) d55_e
READ(8000,*) s5_3_e
READ(8000,*) s56_1_e
READ(8000,*) s56_2_e
READ(8000,*) s56_3_e
READ(8000,*) s3_1_e
READ(8000,*) s3_2_e
READ(8000,*) s3_3_e
READ(8000,*) tencoef22_e
CLOSE(unit=8000)

!---------------------------------------------------------------------------
!(b) Values x in [0,1] are normalized as x*100 so that they lie in [0,100]
!---------------------------------------------------------------------------
p_vector(1)  = DLOG(alpha_1_e/(1.0_dp - alpha_1_e)) 
p_vector(2)  = DLOG(alpha_2_e/(1.0_dp - alpha_2_e))  
p_vector(3)  = DLOG(alpha_3_e/(1.0_dp - alpha_3_e))                                     
p_vector(4)  = DLOG(xi_1_e/(1.0_dp - xi_1_e))
p_vector(5)  = DLOG(xi_2_e/(1.0_dp - xi_2_e))   
p_vector(6)  = DLOG(xi_3_e/(1.0_dp - xi_3_e)) 
p_vector(7)  = gammapast11_0_e
p_vector(8)  = DLOG(p_1_e/(1.0_dp - p_1_e))   
p_vector(9)  = gammap_01_e
p_vector(10) = gammap_02_e
p_vector(11) = gammap_03_e
p_vector(12) = s_1_e
p_vector(13) = s_2_e 
p_vector(14) = DLOG(sigma1_u_e)
p_vector(15) = s1_5_e
p_vector(16) = s_3_e
p_vector(17) = s_002_e
p_vector(18) = s_003_e
p_vector(19) = DLOG(sigma2_u_e)
p_vector(20) = DLOG(sigma3_u_e)
p_vector(21) = w_74_e
p_vector(22) = w_75_e
p_vector(23) = w_76_e
p_vector(24) = w_77_e
p_vector(25) = w_78_e
p_vector(26) = w_79_e
p_vector(27) = gammap_04_e
p_vector(28) = DLOG(d9_e)
p_vector(29) = DLOG(d10_e)
p_vector(30) = DLOG(sigma32_u_e)
p_vector(31) = DLOG(sigma33_u_e)
p_vector(32) = DLOG(sigma34_u_e)
p_vector(33) = DLOG(d5_e)
p_vector(34) = DLOG(phk2_e/(1.0_dp - phk2_e))  
p_vector(35) = DLOG(d12_e)
p_vector(36) = DLOG(xi_32_e/(1.0_dp - xi_32_e)) 
p_vector(37) = DLOG(xi_14_e/(1.0_dp - xi_14_e)) 
p_vector(38) = DLOG(p_3_e/(1.0_dp - p_3_e)) 
p_vector(39) = DLOG(xi_34_e/(1.0_dp - xi_34_e)) 
p_vector(40) = DLOG(d_e)
p_vector(41) = DLOG(d16_e)
p_vector(42) = DLOG(xi_21_e/(1.0_dp - xi_21_e))
p_vector(43) = d1_1_e
p_vector(44) = d1_2_e
p_vector(45) = d1_3_e
p_vector(46) = DLOG(cpen_e/(1.0_dp - cpen_e)) 
p_vector(47) = DLOG(gammapast11_1_e)
p_vector(48) = DLOG(xi_25_e/(1.0_dp - xi_25_e)) 
p_vector(49) = gammapast33_0_e
p_vector(50) = DLOG(d2_e)
p_vector(51) = DLOG(d18_e)
p_vector(52) = DLOG(d19_e)
p_vector(53) = DLOG(d20_e)
p_vector(54) = DLOG(d21_e)
p_vector(55) = DLOG(xi_33_e/(1.0_dp - xi_33_e)) 
p_vector(56) = DLOG(d22_e)
p_vector(57) = DLOG(d23_e)
p_vector(58) = DLOG(d24_e)
p_vector(59) = DLOG(d25_e)
p_vector(60) = DLOG(beta_1_e/(1.0_dp - beta_1_e)) 
p_vector(61) = DLOG(beta_2_e/(1.0_dp - beta_2_e)) 
p_vector(62) = DLOG(beta_3_e/(1.0_dp - beta_3_e))
!p_vector(62) = beta_3_e 
p_vector(63) = DLOG(gammapast33_1_e)
p_vector(64) = DLOG(d26_e)
p_vector(65) = DLOG(d27_e)
p_vector(66) = DLOG(d28_e)
p_vector(67) = d52_e
p_vector(68) = DLOG(d29_e)
p_vector(69) = DLOG(d30_e)
p_vector(70) = DLOG(d31_e)
p_vector(71) = DLOG(delta_e/(1.0_dp - delta_e)) 
p_vector(72) = DLOG(d32_e)
p_vector(73) = DLOG(plus_e/(1.0_dp - plus_e)) 
p_vector(74) = d5j2_e
p_vector(75) = d5j3_e
p_vector(76) = d53_e
p_vector(77) = d533_e
p_vector(78) = DLOG(d523_e)
p_vector(79) = wd_e
p_vector(80) = DLOG(wd3_e)
p_vector(81) = DLOG(wd2_e)
p_vector(82) = DLOG(tencoef_e)
p_vector(83) = DLOG(tencoef1_e)
p_vector(84) = DLOG(tencoef2_e)
p_vector(85) = DLOG(sigma22_u_e)
p_vector(86) = DLOG(sigma23_u_e)
p_vector(87) = DLOG(sigma24_u_e)
p_vector(88) = DLOG(wd52_e)
p_vector(89) = s_0032_e
p_vector(90) = s_0033_e
p_vector(91) = s_0034_e
p_vector(92) = s_0022_e
p_vector(93) = s_0023_e
p_vector(94) = s_0024_e
p_vector(95) = s1_52_e
p_vector(96) = s1_53_e
p_vector(97) = s1_54_e
p_vector(98) = d5t2_e 
p_vector(99) = d5t3_e 
p_vector(100)= d5t4_e
p_vector(101) = DLOG(sigma12_u_e)
p_vector(102) = DLOG(sigma13_u_e)
p_vector(103) = DLOG(sigma14_u_e)
p_vector(104) = wpen_e
p_vector(105) = DLOG(wvar_e)
p_vector(106) = DLOG(v74_e)
p_vector(107) = DLOG(v75_e)
p_vector(108) = DLOG(v76_e)
p_vector(109) = DLOG(v77_e)
p_vector(110) = DLOG(v78_e)
p_vector(111) = DLOG(v79_e)
p_vector(112) = wpen1_e
p_vector(113) = por_e
p_vector(114) = prop_e
p_vector(115) = DLOG(tencoef3_e)
p_vector(116) = tencoef31_e
p_vector(117) = d5j22_e
p_vector(118) = d5j23_e
p_vector(119) = d5j24_e
p_vector(120) = w2_74_e
p_vector(121) = w2_75_e
p_vector(122) = w2_76_e
p_vector(123) = w2_77_e
p_vector(124) = w2_78_e
p_vector(125) = w2_79_e
p_vector(126) = w3_74_e
p_vector(127) = w3_75_e
p_vector(128) = w3_76_e
p_vector(129) = w3_77_e
p_vector(130) = w3_78_e
p_vector(131) = w3_79_e
p_vector(132) = s5_1_e 
p_vector(133) = s5_2_e
p_vector(134) = DLOG(s1_55_e)
p_vector(135) = d55_e 
p_vector(136) = s5_3_e
p_vector(137) = s56_1_e
p_vector(138) = s56_2_e
p_vector(139) = s56_3_e
p_vector(140) = s3_1_e
p_vector(141) = s3_2_e
p_vector(142) = s3_3_e
p_vector(143) = tencoef22_e


!-------------------------------------------------------------
!(c) Indicator function for the parameters actually estimated
!-------------------------------------------------------------
est_ind(1)  = 1 !alpha_1   (prob. high output for high ab. at L1)                       
est_ind(2)  = 1 !alpha_2   (prob. high output for high ab. at L2)                       
est_ind(3)  = 1 !alpha_3   (prob. high output for high ab. at L3)                       
est_ind(4)  = 1 !xi_1      (sep. rate at L1)                                   
est_ind(5)  = 1 !xi_2      (sep. rate at L2)                                   
est_ind(6)  = 1 !xi_3      (sep. rate at L3)                                 
est_ind(7)  = 1 !gammapast11_0                                                          
est_ind(8)  = 1 !p_1       (binomial prop. to generate type prob.)                     
est_ind(9)  = 1 !gammap_01 (parameter for first type prior)                          
est_ind(10) = 1 !gammap_02 (parameter for second type prior)                         
est_ind(11) = 1 !gammap_03 (parameter for third type prior)                          
est_ind(12) = 1 !s_1       (coeff. on age   in wage reg. at L1, L2)                      
est_ind(13) = 1 !s_2       (coeff. on age^2 in wage reg. at L1, L2)                      
est_ind(14) = 1 !sigma1_u  (type-1 wage st.dev at L1)                                   
est_ind(15) = 1 !s1_5      (type-1 const. in wage reg. at L1)              
est_ind(16) = 1 !s_3       (coeff. on edu in wage regression)                            
est_ind(17) = 1 !s_002     (const in wage regression at L2 for type 1)                  
est_ind(18) = 1 !s_003     (const in wage regression at L3 for type 1)                  
est_ind(19) = 1 !sigma2_u  (type-1 wage st.dev at L2)                                   
est_ind(20) = 1 !sigma3_u  (type-1 wage st.dev at L3)                                   
est_ind(21) = 0 !**N0: w_74 = w_75                                                     
est_ind(22) = 1 !w_75      (wage dummy at L1 in 75)                                      
est_ind(23) = 1 !w_76      (wage dummy at L1 in 76)                                     
est_ind(24) = 1 !w_77      (wage dummy at L1 in 77)                                     
est_ind(25) = 1 !w_78      (wage dummy at L1 in 78)                                     
est_ind(26) = 1 !w_79      (wage dummy at L1 in 79)                                     
est_ind(27) = 1 !gammap_04 (parameter for fourth type prior)                         
est_ind(28) = 0 !**NO: d9                         
est_ind(29) = 1 !d10                                               
est_ind(30) = 0 !**NO: sigma32_u (type-2 wage st.dev at L3)                                  
est_ind(31) = 0 !**NO: sigma33_u (type-3 wage st.dev at L3)                                   
est_ind(32) = 0 !**NO: sigma34_u (type-4 wage st.dev at L3)                                   
est_ind(33) = 1 !d5                                            
est_ind(34) = 1 !phk2      (type-2 proprotion)                                          
est_ind(35) = 1 !d12                                                        
est_ind(36) = 0 !**NO: xi_32 = 0                               
est_ind(37) = 1 !xi_14     (sep. rate at L2 and L3)                               
est_ind(38) = 1 !p_3       (type-3 proportion)                                          
est_ind(39) = 1 !xi_34     (sep. rate at at L2 and L3)                           
est_ind(40) = 1 !d                                           
est_ind(41) = 1 !d16                              
est_ind(42) = 1 !xi_21     (sep. rate at L1)                              
est_ind(43) = 1 !d1_1      (coeff. on tenure in class. error for rating at L1)    
est_ind(44) = 1 !d1_2      (coeff. on tenure in class. error for rating at L2)          
est_ind(45) = 0 !**NO: d1_3              
est_ind(46) = 1 !cpen      (sep. rate at L2 and L3)                                
est_ind(47) = 0 !**NO: gammapast11_1 (rating error: d_1(1)=d_1(2))     
est_ind(48) = 1 !xi_25     (sep. rate at  L2 and at L3)                 
est_ind(49) = 0 !**NO: gammapast33_0 (rating error parameter: d_0(2))                         
est_ind(50) = 0 !**NO: d2                                                               
est_ind(51) = 0 !**NO: d18                                        
est_ind(52) = 0 !**NO: d19      
est_ind(53) = 1 !d20       
est_ind(54) = 1 !d21     
est_ind(55) = 0 !**NO: xi_33     
est_ind(56) = 0 !**NO: d22       
est_ind(57) = 1 !d23       
est_ind(58) = 1 !d24    
est_ind(59) = 0 !**NO: d25=0      
est_ind(60) = 1 !beta_1    (prob. high output for high ab. at L1)                      
est_ind(61) = 1 !beta_2    (prob. high output for high ab. at L2)                      
est_ind(62) = 1 !beta_3    (prob. high output for high ab. at L3)                 
est_ind(63) = 0 !**NO: gammapast33_1 (rating error parameter: d_1(2))                        
est_ind(64) = 0 !**NO: d26=0      
est_ind(65) = 1 !d27      
est_ind(66) = 1 !d28      
est_ind(67) = 0 !**NO: d52= d523 
est_ind(68) = 0 !**NO: d29=0                                          
est_ind(69) = 1 !d30      
est_ind(70) = 1 !d31      
est_ind(71) = 0 !**NO: delta (discount factor)                                          
est_ind(72) = 1 !d32             
est_ind(73) = 1 !plus      (sep. rate at L1 and L2)                              
!(LEVEL PARAMETERS UNTIL HERE)
est_ind(74) = 0 !**NO: d5j2=0      
est_ind(75) = 0 !**NO: d5j3=0      
est_ind(76) = 0 !**NO: d53 =0  
est_ind(77) = 0 !**NO: d533 =0 
est_ind(78) = 0 !**NO: d523                                  
est_ind(79) = 0 !**NO: wd=0          
est_ind(80) = 0 !**NO: wd3      
est_ind(81) = 0 !**NO: wd2     
est_ind(82) = 0 !**NO: tencoef = tencoef1      
est_ind(83) = 1 !tencoef1    
est_ind(84) = 0 !**NO: tencoef2      
est_ind(85) = 1 !sigma22_u     (type-2 wage st.dev at L2)                               
est_ind(86) = 1 !sigma23_u     (type-3 wage st.dev at L2)                               
est_ind(87) = 1 !sigma24_u     (type-4 wage st.dev at L2)                               
est_ind(88) = 0 !**NO: wd52      
est_ind(89) = 1 !s_0032        (const in wage regression at L3 for type 2)              
est_ind(90) = 1 !s_0033        (const in wage regression at L3 for type 3)             
est_ind(91) = 1 !s_0034        (const in wage regression at L3 for type 4)              
est_ind(92) = 1 !s_0022        (const in wage regression at L2 for type 2)             
est_ind(93) = 1 !s_0023        (const in wage regression at L2 for type 3)                  
est_ind(94) = 1 !s_0024        (const in wage regression at L2 for type 4)              
est_ind(95) = 1 !s1_52         (type-2 const. in wage reg. at L1 in t=1,2,3,4)          
est_ind(96) = 1 !s1_53         (type-3 const. in wage reg. at L1 in t=1,2,3,4)          
est_ind(97) = 1 !s1_54         (type-4 const. in wage reg. at L1 in t=1,2,3,4)         
est_ind(98) = 1 !d5t2          (type-2 slope on p at L1)                               
est_ind(99) = 1 !d5t3          (type-3 slope on p at L1)                                
est_ind(100)= 1 !d5t4          (type-4 slope on p at L1)                                
est_ind(101)= 1 !sigma12_u     (type-2 wage st.dev at L1)                               
est_ind(102)= 1 !sigma13_u     (type-3 wage st.dev at L1)                               
est_ind(103)= 1 !sigma14_u     (type-4 wage st.dev at L1)                               
est_ind(104)= 0 !**NO: wpen=0    
est_ind(105)= 0 !**NO: wvar=0  
est_ind(106)= 0 !**NO: v74=1   
est_ind(107)= 0 !**NO: v75=1     
est_ind(108)= 0 !**NO: v76=1     
est_ind(109)= 0 !**NO: v77=1      
est_ind(110)= 0 !**NO: v78=1  
est_ind(111)= 0 !**NO: v79=1      
est_ind(112)= 0 !**NO: wpen1=0 
est_ind(113)= 0 !**NO: por=0   
est_ind(114)= 0 !**NO: prop          
est_ind(115)= 0 !**NO: tencoef3      
est_ind(116)= 0 !**NO: tencoef31    
est_ind(117)= 0 !**NO: d5j22 =d5j3 
est_ind(118)= 0 !**NO: d5j23 =d5j3    
est_ind(119)= 0 !**NO: d5j24 =d5j3
est_ind(120)= 0 !**NO: w2_74 = w2_75 (wage dummy at L2 in 74)                          
est_ind(121)= 0 !**NO: w2_75         (wage dummy at L2 in 75)                                 
est_ind(122)= 0 !**NO: w2_76         (wage dummy at L2 in 76)                                 
est_ind(123)= 0 !**NO: w2_77         (wage dummy at L2 in 77)                                
est_ind(124)= 0 !**NO: w2_78         (wage dummy at L2 in 78)                                
est_ind(125)= 0 !**NO: w2_79         (wage dummy at L2 in 79)                                 
est_ind(126)= 0 !**NO: w3_74 = w3_75 (wage dummy at L3 in 74)                            
est_ind(127)= 0 !**NO: w3_75         (wage dummy at L3 in 75)                                
est_ind(128)= 0 !**NO: w3_76         (wage dummy at L3 in 76)                                 
est_ind(129)= 0 !**NO: w3_77         (wage dummy at L3 in 77)                                  
est_ind(130)= 0 !**NO: w3_78         (wage dummy at L3 in 78)                                  
est_ind(131)= 0 !**NO: w3_79         (wage dummy at L3 in 79)                                 
est_ind(132)= 0 !**NO: s5_1=0 
est_ind(133)= 0 !**NO: s5_2=0 
est_ind(134)= 0 !**NO: s1_55   
est_ind(135)= 0 !**NO: d55=0  
est_ind(136)= 0 !**NO: s5_3=0  
est_ind(137)= 0 !**NO: s56_1=0 
est_ind(138)= 0 !**NO: s56_2=0 
est_ind(139)= 0 !**NO: s56_3=0 
est_ind(140)= 1 !s3_1          (coeff. on age    in wage reg. at L3)                     
est_ind(141)= 1 !s3_2          (coeff. on age**2 in wage reg. at L3)                    
est_ind(142)= 1 !s3_3          (coeff. on edu    in wage reg. at L3)                     
est_ind(143)= 0 !**NO: tencoef22     


!********************************************
!(3) Number of parameters to be estimated
!********************************************
num_param = SUM(est_ind)


!****************************************************
!(4) File reporting which parameters are estimated
!****************************************************
OPEN(unit=7777, file='estimated_parameters.txt', status='replace')
WRITE(7777,*)  "Number of Parameters to be Estimated =",num_param
WRITE(7777,*) "Total Estimated Parameters:"

IF (est_ind(1)==1) THEN
  WRITE(7777,*) "(1) Alpha_1"
END IF

IF (est_ind(2)==1) THEN
  WRITE(7777,*) "(2) Alpha_2"
END IF

IF (est_ind(3)==1) THEN
  WRITE(7777,*) "(3) Alpha_3"
END IF

IF (est_ind(4)==1) THEN
  WRITE(7777,*) "(4) xi_1"
END IF

IF (est_ind(5)==1) THEN
  WRITE(7777,*) "(5) xi_2"
END IF

IF (est_ind(6)==1) THEN
  WRITE(7777,*) "(6) xi_3"
END IF

IF (est_ind(7)==1) THEN
  WRITE(7777,*) "(7) gammapast11_0"
END IF

IF (est_ind(8)==1) THEN
  WRITE(7777,*) "(8) p_1"
END IF

IF (est_ind(9)==1) THEN
  WRITE(7777,*) "(9)  gammap_01"
END IF

IF (est_ind(10)==1) THEN
  WRITE(7777,*) "(10)  gammap_02"
END IF

IF (est_ind(11)==1) THEN
  WRITE(7777,*) "(11)  gammap_03"
END IF

 IF (est_ind(12)==1) THEN
  WRITE(7777,*) "(12) s_1"
END IF

IF (est_ind(13)==1) THEN
  WRITE(7777,*) "(13) s_2"
END IF

IF (est_ind(14)==1) THEN
  WRITE(7777,*) "(14) sigma1_u"
END IF

IF (est_ind(15)==1) THEN
  WRITE(7777,*) "(15) s1_5"
END IF

IF (est_ind(16)==1) THEN
  WRITE(7777,*) "(16) s_3"
END IF

IF (est_ind(17)==1) THEN
  WRITE(7777,*) "(17) s_002"
END IF

IF (est_ind(18)==1) THEN
  WRITE(7777,*) "(18) s_003"
END IF

IF (est_ind(19)==1) THEN
  WRITE(7777,*) "(19) sigma2_u"
END IF

IF (est_ind(20)==1) THEN
  WRITE(7777,*) "(20) sigma3_u"
END IF

IF (est_ind(21)==1) THEN
  WRITE(7777,*) "(21) w_74"
END IF

IF (est_ind(22)==1) THEN
  WRITE(7777,*) "(22) w_75"
END IF

IF (est_ind(23)==1) THEN
  WRITE(7777,*) "(23) w_76"
END IF

IF (est_ind(24)==1) THEN
  WRITE(7777,*) "(24) w_77"
END IF

IF (est_ind(25)==1) THEN
  WRITE(7777,*) "(25) w_78"
END IF

IF (est_ind(26)==1) THEN
  WRITE(7777,*) "(26) w_79"
END IF

IF (est_ind(27)==1) THEN
  WRITE(7777,*) "(27) gammap_04"
END IF

IF (est_ind(28)==1) THEN
  WRITE(7777,*) "(28) d9"
END IF

IF (est_ind(29)==1) THEN
  WRITE(7777,*) "(29) d10"
END IF

IF (est_ind(30)==1) THEN
  WRITE(7777,*) "(30) sigma32_u"
END IF

IF (est_ind(31)==1) THEN
  WRITE(7777,*) "(31) sigma33_u"
END IF

IF (est_ind(32)==1) THEN
  WRITE(7777,*) "(32) sigma34_u"
END IF

IF (est_ind(33)==1) THEN
  WRITE(7777,*) "(33) d5"
END IF

IF (est_ind(34)==1) THEN
  WRITE(7777,*) "(34) phk2"
END IF

IF (est_ind(35)==1) THEN
  WRITE(7777,*) "(35) d12"
END IF

IF (est_ind(36)==1) THEN
  WRITE(7777,*) "(36) xi_32"
END IF

IF (est_ind(37)==1) THEN
  WRITE(7777,*) "(37) xi_14"
END IF

IF (est_ind(38)==1) THEN
  WRITE(7777,*) "(38) p_3"
END IF

IF (est_ind(39)==1) THEN
  WRITE(7777,*) "(39) xi_34"
END IF

IF (est_ind(40)==1) THEN
  WRITE(7777,*) "(40) d"
END IF

IF (est_ind(41)==1) THEN
  WRITE(7777,*) "(41) d16"
END IF

IF (est_ind(42)==1) THEN
  WRITE(7777,*) "(42) xi_21"
END IF

IF (est_ind(43)==1) THEN
  WRITE(7777,*) "(43) d1_1"
END IF

IF (est_ind(44)==1) THEN
  WRITE(7777,*) "(44) d1_2"
END IF

IF (est_ind(45)==1) THEN
  WRITE(7777,*) "(45) d1_3"
END IF

IF (est_ind(46)==1) THEN
  WRITE(7777,*) "(46) cpen"
END IF
IF (est_ind(47)==1) THEN
  WRITE(7777,*) "(47)  gammapast11_1"
END IF
IF (est_ind(48)==1) THEN
  WRITE(7777,*) "(48) xi_25"
END IF
IF (est_ind(49)==1) THEN
  WRITE(7777,*) "(49)  gammapast33_0"
END IF
IF (est_ind(50)==1) THEN
  WRITE(7777,*) "(50) d2"
END IF
IF (est_ind(51)==1) THEN
  WRITE(7777,*) "(51) d18"
END IF
IF (est_ind(52)==1) THEN
  WRITE(7777,*) "(52) d19"
END IF
IF (est_ind(53)==1) THEN
  WRITE(7777,*) "(53) d20"
END IF
IF (est_ind(54)==1) THEN
  WRITE(7777,*) "(54) d21"
END IF
IF (est_ind(55)==1) THEN
  WRITE(7777,*) "(55) xi_33"
END IF
IF (est_ind(56)==1) THEN
  WRITE(7777,*) "(56) d22"
END IF
IF (est_ind(57)==1) THEN
  WRITE(7777,*) "(57) d23"
END IF
IF (est_ind(58)==1) THEN
  WRITE(7777,*) "(58) d24"
END IF
IF (est_ind(59)==1) THEN
  WRITE(7777,*) "(59) d25"
END IF
IF (est_ind(60)==1) THEN
  WRITE(7777,*) "(60) beta_1"
END IF
IF (est_ind(61)==1) THEN
  WRITE(7777,*) "(61) beta_2"
END IF
IF (est_ind(62)==1) THEN
  WRITE(7777,*) "(62) beta_3"
END IF
IF (est_ind(63)==1) THEN
  WRITE(7777,*) "(63)  gammapast33_1"
END IF
IF (est_ind(64)==1) THEN
  WRITE(7777,*) "(64) d26"
END IF
IF (est_ind(65)==1) THEN
  WRITE(7777,*) "(65) d27"
END IF
IF (est_ind(66)==1) THEN
  WRITE(7777,*) "(66) d28"
END IF
IF (est_ind(67)==1) THEN
  WRITE(7777,*) "(67) d52"
END IF
IF (est_ind(68)==1) THEN
  WRITE(7777,*) "(68) d29"
END IF
IF (est_ind(69)==1) THEN
  WRITE(7777,*) "(69) d30"
END IF
IF (est_ind(70)==1) THEN
  WRITE(7777,*) "(70) d31"
END IF
IF (est_ind(71)==1) THEN
  WRITE(7777,*) "(71) delta"
END IF
IF (est_ind(72)==1) THEN
  WRITE(7777,*) "(72) d32"
END IF
IF (est_ind(73)==1) THEN
  WRITE(7777,*) "(73) plus"
END IF
IF (est_ind(74)==1) THEN
  WRITE(7777,*) "(74) d5j2"
END IF
IF (est_ind(75)==1) THEN
  WRITE(7777,*) "(75) d5j3"
END IF
IF (est_ind(76)==1) THEN
  WRITE(7777,*) "(76) d53"
END IF
IF (est_ind(77)==1) THEN
  WRITE(7777,*) "(77) d533"
END IF
IF (est_ind(78)==1) THEN
  WRITE(7777,*) "(78) d523"
END IF
IF (est_ind(79)==1) THEN
  WRITE(7777,*) "(79) wd"
END IF
IF (est_ind(80)==1) THEN
  WRITE(7777,*) "(80) wd3"
END IF
IF (est_ind(81)==1) THEN
  WRITE(7777,*) "(81) wd2"
END IF
IF (est_ind(82)==1) THEN
  WRITE(7777,*) "(82) tencoef"
END IF
IF (est_ind(83)==1) THEN
  WRITE(7777,*) "(83) tencoef1"
END IF
IF (est_ind(84)==1) THEN
  WRITE(7777,*) "(84) tencoef2"
END IF
IF (est_ind(85)==1) THEN
  WRITE(7777,*) "(85) sigma22_u"
END IF
IF (est_ind(86)==1) THEN
  WRITE(7777,*) "(86) sigma23_u"
END IF
IF (est_ind(87)==1) THEN
  WRITE(7777,*) "(87) sigma24_u"
END IF
IF (est_ind(88)==1) THEN
  WRITE(7777,*) "(88) wd52"
END IF
IF (est_ind(89)==1) THEN
  WRITE(7777,*) "(89) s_0032"
END IF
IF (est_ind(90)==1) THEN
  WRITE(7777,*) "(90) s_0033"
END IF
IF (est_ind(91)==1) THEN
  WRITE(7777,*) "(91) s_0034"
END IF
IF (est_ind(92)==1) THEN
  WRITE(7777,*) "(92) s_0022"
END IF
IF (est_ind(93)==1) THEN
  WRITE(7777,*) "(93) s_0023"
END IF
IF (est_ind(94)==1) THEN
  WRITE(7777,*) "(94) s_0024"
END IF
IF (est_ind(95)==1) THEN
  WRITE(7777,*) "(95) s1_52"
END IF
IF (est_ind(96)==1) THEN
  WRITE(7777,*) "(96) s1_53"
END IF
IF (est_ind(97)==1) THEN
  WRITE(7777,*) "(97) s1_54"
END IF
IF (est_ind(98)==1) THEN
  WRITE(7777,*) "(98) d5t2"
END IF
IF (est_ind(99)==1) THEN
  WRITE(7777,*) "(99) d5t3"
END IF
IF (est_ind(100)==1) THEN
  WRITE(7777,*) "(100) d5t4"
END IF
IF (est_ind(101)==1) THEN
  WRITE(7777,*) "(101) sigma12_u"
END IF
IF (est_ind(102)==1) THEN
  WRITE(7777,*) "(102) sigma13_u"
END IF
IF (est_ind(103)==1) THEN
  WRITE(7777,*) "(103) sigma14_u"
END IF
IF (est_ind(104)==1) THEN
  WRITE(7777,*) "(104) wpen"
END IF
IF (est_ind(105)==1) THEN
  WRITE(7777,*) "(105) wvar"
END IF
IF (est_ind(106)==1) THEN
  WRITE(7777,*) "(106) v74"
END IF

IF (est_ind(107)==1) THEN
  WRITE(7777,*) "(107) v75"
END IF

IF (est_ind(108)==1) THEN
  WRITE(7777,*) "(108) v76"
END IF

IF (est_ind(109)==1) THEN
  WRITE(7777,*) "(109) v77"
END IF

IF (est_ind(110)==1) THEN
  WRITE(7777,*) "(110) v78"
END IF

IF (est_ind(111)==1) THEN
  WRITE(7777,*) "(111) v79"
END IF
IF (est_ind(112)==1) THEN
  WRITE(7777,*) "(112) wpen1"
END IF
IF (est_ind(113)==1) THEN
  WRITE(7777,*) "(113) por"
END IF
IF (est_ind(114)==1) THEN
  WRITE(7777,*) "(114) prop"
END IF
IF (est_ind(115)==1) THEN
  WRITE(7777,*) "(115) tencoef3"
END IF
IF (est_ind(116)==1) THEN
  WRITE(7777,*) "(116) tencoef31"
END IF
IF (est_ind(117)==1) THEN
  WRITE(7777,*) "(117) d5j22"
END IF
IF (est_ind(118)==1) THEN
  WRITE(7777,*) "(118) d5j23"
END IF
IF (est_ind(119)==1) THEN
  WRITE(7777,*) "(119) d5j24"
END IF
IF (est_ind(120)==1) THEN
  WRITE(7777,*) "(120) w2_74"
END IF
IF (est_ind(121)==1) THEN
  WRITE(7777,*) "(121) w2_75"
END IF
IF (est_ind(122)==1) THEN
  WRITE(7777,*) "(122) w2_76"
END IF
IF (est_ind(123)==1) THEN
  WRITE(7777,*) "(123) w2_77"
END IF
IF (est_ind(124)==1) THEN
  WRITE(7777,*) "(124) w2_78"
END IF
IF (est_ind(125)==1) THEN
  WRITE(7777,*) "(125) w2_79"
END IF
IF (est_ind(126)==1) THEN
  WRITE(7777,*) "(126) w3_74"
END IF
IF (est_ind(127)==1) THEN
  WRITE(7777,*) "(127) w3_75"
END IF
IF (est_ind(128)==1) THEN
  WRITE(7777,*) "(128) w3_76"
END IF
IF (est_ind(129)==1) THEN
  WRITE(7777,*) "(129) w3_77"
END IF
IF (est_ind(130)==1) THEN
  WRITE(7777,*) "(130) w3_78"
END IF
IF (est_ind(131)==1) THEN
  WRITE(7777,*) "(131) w3_79"
END IF
IF (est_ind(132)==1) THEN
  WRITE(7777,*) "(132) s5_1"
END IF
IF (est_ind(133)==1) THEN
  WRITE(7777,*) "(133) s5_2"
END IF
IF (est_ind(134)==1) THEN
  WRITE(7777,*) "(134) s1_55"
END IF
IF (est_ind(135)==1) THEN
  WRITE(7777,*) "(135) d55"
END IF
IF (est_ind(136)==1) THEN
  WRITE(7777,*) "(136) s5_3"
END IF

IF (est_ind(137)==1) THEN
  WRITE(7777,*) "(137) s56_1"
END IF

IF (est_ind(138)==1) THEN
  WRITE(7777,*) "(138) s56_2"
END IF
IF (est_ind(139)==1) THEN
  WRITE(7777,*) "(139) s56_3"
END IF
IF (est_ind(140)==1) THEN
  WRITE(7777,*) "(140) s3_1"
END IF
IF (est_ind(141)==1) THEN
  WRITE(7777,*) "(141) s3_2"
END IF
IF (est_ind(142)==1) THEN
  WRITE(7777,*) "(142) s3_3"
END IF
IF (est_ind(143)==1) THEN
  WRITE(7777,*) "(143) tencoef22"
END IF


CLOSE (unit=7777)

END SUBROUTINE param_initialization



!********************************************************************************************************
!(5) Subroutine to assign the initial values to the P parameter matrix (input for the amoeba subroutine)
!********************************************************************************************************
SUBROUTINE matrix_initial(pr)

IMPLICIT NONE

!REMARK: output is the parameter matrix P to be estimated and the correspondent position in the parameter vector p_vector


REAL(dp), DIMENSION(num_param+1,num_param), INTENT(OUT):: pr
INTEGER:: j1,j2
REAL(dp), DIMENSION(tot_param):: lamda
INTEGER, PARAMETER:: comp_sterr = 1
REAL(dp):: s = 0.025_dp



ALLOCATE(pos(num_param))


!----------------------------------------------------------------------------------------------------------------------------------
!(a) Vector of bumps (lamda is the perturbation vector used for constructing the param_matrix used in loglikelihood maximization)
!----------------------------------------------------------------------------------------------------------------------------------
!-----------------------------
!i. For regular computations
!-----------------------------
IF (comp_sterr == 0) THEN
lamda(1) = 0.20_dp
lamda(2) = 0.20_dp
lamda(3) = 0.20_dp
lamda(4) = 0.015_dp  
lamda(5) = 0.015_dp 
lamda(6) = 0.015_dp  
lamda(7) = 0.35_dp  
lamda(8) = 0.25_dp
lamda(9) = 0.15_dp
lamda(10)= 0.40_dp 
lamda(11)= -0.60_dp
lamda(12)= 0.15_dp
lamda(13)= 0.15_dp
lamda(14)= -0.15_dp
lamda(15)= 0.005_dp
lamda(16)= 0.30_dp 
lamda(17)= 0.05_dp  
lamda(18)= 0.05_dp   
lamda(19)= 0.20_dp  
lamda(20)= -0.15_dp
lamda(21)= 0.30_dp
lamda(22)= 0.30_dp
lamda(23)= 0.30_dp
lamda(24)= 0.30_dp
lamda(25)= 0.30_dp
lamda(26)= 0.30_dp
lamda(27)= 0.55_dp
lamda(28)= 0.005_dp
lamda(29)= 0.005_dp
lamda(30)= 0.25_dp
lamda(31)= 0.25_dp
lamda(32)= 0.25_dp
lamda(33)= 0.05_dp
lamda(34)= 0.25_dp
lamda(35)= 0.0015_dp
lamda(36)= 0.60_dp
lamda(37)= 0.015_dp
lamda(38)= 0.35_dp
lamda(39)= 0.015_dp
lamda(40)= 0.015_dp
lamda(41)= 0.015_dp
lamda(42)= 0.015_dp
lamda(43)= 0.15_dp
lamda(44)= 0.15_dp
lamda(45)= 0.20_dp
lamda(46)= 0.015_dp
lamda(47)= 0.35_dp
lamda(48)= 0.015_dp
lamda(49)= 0.30_dp
lamda(50)= 0.015_dp
lamda(51)= 0.20_dp
lamda(52)= 0.15_dp
lamda(53)= 0.0015_dp
lamda(54)= 0.0015_dp
lamda(55)= 0.30_dp
lamda(56)= 0.0015_dp
lamda(57)= 0.0015_dp
lamda(58)= 0.05_dp
lamda(59)= 0.015_dp
lamda(60)= 0.10_dp
lamda(61)= 0.35_dp
lamda(62)= 0.35_dp
lamda(63)= 0.35_dp
lamda(64)= 0.015_dp
lamda(65)= 0.015_dp
lamda(66)= 0.005_dp
lamda(67)= 0.10_dp
lamda(68)= 0.10_dp
lamda(69)= 0.015_dp
lamda(70)= 0.015_dp
lamda(71)= 0.15_dp
lamda(72)= 0.0015_dp
lamda(73)= 0.15_dp
lamda(74)= 0.25_dp
lamda(75)= 0.15_dp
lamda(76)= 0.15_dp
lamda(77)= 0.25_dp
lamda(78)= 1.25_dp
lamda(79)= 1.25_dp
lamda(80)= 1.25_dp
lamda(81)= 1.25_dp
lamda(82)= 0.65_dp
lamda(83)= 0.35_dp
lamda(84)= 0.95_dp
lamda(85)= 0.25_dp
lamda(86)= 0.25_dp
lamda(87)= 0.25_dp
lamda(88)= 0.65_dp
lamda(89)= 0.25_dp
lamda(90)= 0.25_dp
lamda(91)= 0.25_dp
lamda(92)= 0.25_dp
lamda(93)= 0.25_dp
lamda(94)= 0.25_dp
lamda(95)= 0.25_dp
lamda(96)= 0.25_dp
lamda(97)= 0.25_dp
lamda(98)= 0.10_dp
lamda(99)= 0.10_dp
lamda(100)=0.10_dp
lamda(101)=0.10_dp
lamda(102)=0.10_dp
lamda(103)=0.10_dp
lamda(104)=0.20_dp
lamda(105)=0.10_dp
lamda(106)=0.20_dp
lamda(107)=0.20_dp
lamda(108)=0.20_dp
lamda(109)=0.20_dp
lamda(110)=0.20_dp
lamda(111)=0.20_dp
lamda(112)=0.20_dp
lamda(113)=-0.20_dp
lamda(114)= 0.20_dp
lamda(115)= 0.25_dp
lamda(116)= 0.25_dp
lamda(117)= 0.20_dp
lamda(118)= 0.20_dp
lamda(119)= 0.20_dp
lamda(120)= 0.30_dp
lamda(121)= 0.30_dp
lamda(122)= 0.30_dp
lamda(123)= 0.30_dp
lamda(124)= 0.30_dp
lamda(125)= 0.30_dp
lamda(126)= 0.30_dp
lamda(127)= 0.30_dp
lamda(128)= 0.30_dp
lamda(129)= 0.30_dp
lamda(130)= 0.30_dp
lamda(131)= 0.30_dp
lamda(132)= 0.30_dp
lamda(133)= 0.30_dp
lamda(134)= 0.30_dp
lamda(135)= 0.30_dp
lamda(136)= 0.30_dp
lamda(137)= 0.30_dp
lamda(138)= 0.30_dp
lamda(139)= 0.30_dp
lamda(140)= 0.30_dp
lamda(141)= 0.30_dp
lamda(142)= 0.30_dp
lamda(143)= 0.20_dp


!--------------------------------------------------------------------------------------------------------
!ii. For standard error computations (s=0.025 is the perturbation magnitude for computing standard error)
!--------------------------------------------------------------------------------------------------------


ELSE IF (comp_sterr ==1) THEN
DO j1 = 1, tot_param

 lamda(j1)= s

END DO


ELSE

WRITE(3000,*) "Program stopped bc unclear whether or not st errors are to be computed"
STOP

!------------------------------------------
END IF !for computation of standard errors
!------------------------------------------

!-----------------------------------------
!(b) First row of the parameter matrix P
!-----------------------------------------
!NOTE: The vector 'pos' keeps track of the position of the estimated parameter in the ORIGINAL parameter vector
j2 = 1
DO j1 = 1, tot_param
  IF (est_ind(j1)==1) THEN 
        pr(1,j2) = p_vector(j1)
        pos(j2) = j1
        j2 = j2+1
  END IF
END DO




!---------------------------------------------------
!(c) Second to last row of the parameter matrix P
!---------------------------------------------------
DO j1 = 2, num_param+1           !rows of matrix
 DO j2 = 1, num_param            !columns of matrix
   IF(j2==j1-1) THEN
         pr(j1,j2) = pr(1,j2)*(1.0_dp+lamda(pos(j2)))
   ELSE
         pr(j1,j2) = pr(1,j2)
   END IF
 END DO
END DO



!DEALLOCATE(pos)
END SUBROUTINE matrix_initial

END MODULE mod_param_initial




